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M33 X-7 is among the most massive X-Ray binary stellar systems known, hosting a rapidly 
spinning 15.65 M black hole orbiting an underluminous 70 M Main Sequence companion 
in a slightly eccentric 3.45 day orbit 1 ' 2 . Although post-main-sequence mass transfer explains 
the masses and tight orbit 3 , it leaves unexplained the observed X-Ray luminosity, star's un- 
derluminosity, black hole's spin, and eccentricity. A common envelope phase 1 , or rotational 
mixing 4 , could explain the orbit, but the former would lead to a merger and the latter to an 
overluminous companion. A merger would also ensue if mass transfer to the black hole were 
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invoked for its spin-up 5 . Here we report that, if M33 X-7 started as a primary of 85-99 M 
and a secondary of 28-32 M , in a 2.8-3.1 day orbit, its observed properties can be consis- 
tently explained. In this model, the Main Sequence primary transferred part of its envelope 
to the secondary and lost the rest in a wind; it ended its life as a ~16 M He star with a Fe-Ni 
core which collapsed to a black hole (with or without an accompanying supernova). The re- 
lease of binding energy and, possibly, collapse asymmetries "kicked" the nascent black hole 
into an eccentric orbit. Wind accretion explains the X-Ray luminosity, while the black hole 
spin can be natal. 

M33 X-7 has been identified as an evolutionary challenge, given the massive components 
and its tight orbit relative to the large H-rich black hole (BH) companion. Four paths have been 
proposed to explain M33 X-7's formation, but none of them has addressed nor can simultaneously 
explain all its observed properties (Tables 1 and 2). We performed detailed binary evolution cal- 
culations to explore possible evolutionary tracks. Given the spatial metallicity gradient of the M33 
galaxy 6 , we assume a metallicity 50% of the solar value for all our models. 

To illustrate clearly M33 X-7's evolutionary history, we show the results for one of the suc- 
cessful evolutionary sequences in Figure 1. The progenitor comprises a primary of ~97M (BH 
progenitor) and a secondary of ~32 M (BH-companion progenitor) in an orbit of ~2.9 days. Dur- 
ing the first ~1.8 Myr the evolution is driven by mass loss via stellar winds, causing a decrease 
of the gravitational attraction between the components and expansion of the orbit to ~3. 25 days. 
The more massive primary evolves faster than the secondary, growing in size to accommodate the 
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energy produced by fusing H into He at its center. Eventually, while still on its Main Sequence, it 
expands and begins mass transfer (MT) onto the secondary when entering the sphere of influence 
of the secondary's gravitational field (through Roche-lobe overflow). This stronger mode of mass 
loss brings the primary out of thermal equilibrium; in response the star shrinks, recovering its ther- 
mal equilibrium while always maintaining hydrostatic equilibrium, and hence dynamical stability. 
During the first few tens of thousands of years of MT, the orbital period decreases because the more 
massive primary is transferring mass to the less massive secondary. When the secondary accretes 
enough matter to become the more massive component, the orbit starts expanding 7 . The primary 
transfers most of its H-rich envelope and becomes a Wolf-Rayet star, and the strong Wolf-Rayet 
wind (~2 to 3-1CT 5 M yr~ 1 ) removes much of the remaining envelope, eventually interrupting the 
MT. During the ~99,000 years of conservative MT, the original 32 M secondary becomes a mas- 
sive ~69 M O-type star, while the primary becomes a ~51 M Wolf-Rayet. Once the Wolf-Rayet 
wind sets in and the MT is interrupted, the wind blows away the remaining primary's envelope to 
expose the ~25 M He core. This mass loss drives further orbital expansion until the end of the 
primary's Main Sequence and throughout the core He burning phase. At the same time, the now 
more massive secondary is losing mass via its own O-star wind at a lower rate (~ 1CT 6 M yr _1 ). 
At this time the orbit of the binary is circular and the spin period of each star is expected to be syn- 
chronized with the orbital period. The synchronization is due to exchange of angular momentum 
between the stars and their orbit caused by tidal interaction. 

The final stages of the primary's life during and beyond carbon burning are too short (~60 yr 
for an initially ~25 M He star) 8 to significantly change the stellar and orbital parameters. At the 
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end of the primary's life, after ~3.7 Myr, M33 X-7 comprises a ~16 M evolved Wolf-Rayet star 
with an Fe-Ni core, and a ~64.5 M O-star companion in a ~3.5 day orbit. Unable to support itself 
through further nuclear fusion, this massive He-rich star collapses into a BH, and a small fraction of 
the rest mass energy (10%) is released as the BH's gravitational energy. Additionally, asymmetries 
in the collapse and associated neutrino emission may impart an instantaneous linear momentum 
recoil (kick) to the newly born BH, even without any baryonic mass ejection at collapse. Both 
these effects modify the orbital configuration inducing an eccentricity, and slightly decreasing the 
orbital separation to ~3.4days. In fact, while the release of binding energy leads to an increase of 
the orbital separation, the kick imparted to the BH acts to shrink it. For the remaining ~0.2 Myr, 
the post-BH-formation binary evolution is driven by mass loss via the secondary's stellar wind, 
causing the orbit to further expand to the currently observed value. The fraction of this stellar wind 
attracted and accreted by the BH is too small to significantly influence the orbital evolution, but it 
is adequate to explain the X-Ray luminosity observed. At the present time, after ~3.9Myr, M33 
X-7 comprises a BH of ~14.4M and an underluminous O-star of ~64M orbiting around their 
common center of mass in a slightly eccentric ~3.45 day orbit. 

Our model is consistent with a natal nature of the BH's observationally inferred high spin 9 . In 
fact, although it has been suggested that such a high spin is the result of a MT from the companion 
star to the BH through Roche-lobe overflow 5 , given how more massive is the BH companion com- 
pared to the BH, such a phase could not have occurred: it would have been dynamically unstable 
and rapidly evolve into a merger of the binary components. Wind accretion is too weak to spin-up 
the BH to the current value if it were born spinning much slower. In our model, when the BH pro- 
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genitor leaves the Main Sequence, the spins of the stars are expected to be synchronized with the 
orbit, and, assuming solid-body rotation on the Main Sequence, the inner parts of the primary's core 
carry enough angular momentum to explain the currently observed BH spin (the inner 15.5 M of 
the core carry ~5.2x 10 51 grcm 2 s _1 of angular momentum, while ~1.8x 10 51 grcm 2 s _1 is needed 
to explain the currently observed BH's spin). At the end of the Main Sequence, the inner core of 
the star is expected to rotationally decouple from the outer envelope and approximately retain the 
angular momentum of the central layers (ref. 10 , Table 4). 

We explore binary evolutionary sequences with different combinations of initial masses and 
orbital periods. We select as "successful" sequences those that eventually match all observed prop- 
erties within la errors. All these sequences follow a path qualitatively very similar to the specific 
example described in detail here. The progenitors are constrained to host 96-99 M primaries, and 
32 M (within 1 M uncertainty) secondaries in orbits with initial periods of 2.8-2.9 days. The ap- 
parent puzzling underluminosity of the BH companion is due to two factors: (i) the orientation of 
the system with respect to our line of sight and associated projection effects reduce the star's mea- 
sured luminosity (accounting for ~87% of the underluminosity); (ii) the secondary was not born 
as a ~63-65 M star, but instead accreted much of its mass from the BH progenitor (accounting 
for ~13% of the underluminosity) (see Supplementary Information). 

We note that the distance to M33 is more uncertain than 840±20 kpc, and some of the 
observed system properties vary if a different distance to the system is considered (see Table 1). 
Various studies in the literature have reported it in the range 750-1017 kpc (see Supplementary 
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Information). Considering this full range the progenitors are constrained to host primaries between 
85-99 M , secondaries between 28-32 M , and initial orbital periods between 2.8-3.1 days (see 
Figure 2). 

The different phases we described for M33 X-7's evolutionary past have been observed in 
a variety of other binary systems. For example, LMC R136-38 hosts two non-interacting O-stars 
of ~57M and ~23M orbiting around each other every ~3.4days n . LMC-SC1-105 comprises 
a ~31 M O-star that is transferring mass to its ~13M companion 12 . In the system WR46, a 
~51M Wolf-Rayet star orbits a ~6OM O-star companion every ~6days. In particular, the 
various discoveries of Wolf-Rayet stars with O-star companions 13 validate our theoretical model, 
as they resemble the configuration of M33 X-7 less than 2 Myr ago. Several of these binaries have 
been observed with orbital periods of a few days, and massive components (up to ~8OM ) 14,15 . 
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Parameter 


Value 


Parameter 


Value 


M BH (M ) 


15.65 ± 1.45 


P (days) 


3.45301 ± 0.00002 


M 2 (M ) 


70.0 ± 6.9 


e 


0.0185 ± 0.0077 


Spectral Type 


07 III to 08 III 


<(") 


74.6 ± 1.0 


Ted (K) 


35000 ± 1000 


L X (10 38 erg s- 1 ) 


0.13 to 2.49 


Log(L 2 /L ) 


5.72 ± 0.07 


a* 


0.84 ± 0.05 
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Table 1 : Observed parameters for M33 X-7 for a reference distance of 840±20 kpc as 
adopted by the discovery team 1 . The BH mass (M BH ), its companion mass (M 2 ), Spec- 
tral Type, effective temperature (T eff ), and luminosity (L 2 ), the orbital eccentricity (e), and 
inclination (i) are listed as reported by ref. 1 . The orbital period (P) has been measured 
by ref. 2 . The dimensionless spin parameter of the BH (a*) has been determined by ref. 
9 based on the X-ray continuum fitting method 16-18 . The X-ray luminosity (L x ) is derived 
from observations reported in refs. 1 9> 2 °> 2 > 21 > 1 >9. To account for variations in the X-Ray 
flux over different observations, we consider the lowest and highest reported values, after 
we rescale each L x to a M33 distance of 840±20 kpc. If the full distance range of 750- 
101 7 kpc is adopted, using the ELC code of Orosz and Hauschildt 22 we calculate that 
the masses are between 55-1 03 M and 1 3.5-20 M for the star and the BH, respectively, 
and the inclination is between 77°-71°. The logarithmic luminosity in solar units is then 
between 5.62-5. 89 1 . For each distance, the mass of the star can be derived from M 2 = 
-75.94 + 0.17 • d, and the inclination from % = 93.69 - 0.02 • d, where d is in kpc, M 2 in 
M and i in degrees. For each M 2 , the corresponding BH's mass in solar units can be 
calculated from M BH = 6.19 + 0.13 M 2 , and the BH's spin from a* = 0.31 + 0.049 M BH - 
0.001 Mg H . The rescaled X-Ray luminosity ranges from ~1 x 1 37 to ~3x 1 38 erg sec 1 . 
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Table 2: Models suggested for M33 X-7's formation. The various parameters are de- 
scribed in Table 1 . The symbol V" means a parameter has been addressed and ex- 
plained, "x" a parameter has been addressed but not explained, and "O" a parameter 
has not been addressed at all. To explain the tight orbit, ref. 1 suggests that the pro- 
genitors underwent a common envelope (CE) phase during which the primary expanded 
to the point of engulfing its companion. Such a phase is known to lead to tight systems 
because energy can be transferred from the orbit to the CE, leading to a reduction of the 
binary separation, and ejection of the envelope 23 . To form the observed BH, ref. 1 requires 
that the CE begins after He core burning in the primary is completed. However, a CE in 
M33 X-7's case would likely evolve into a merger, because massive-star envelopes are 
tightly bound 24 . Furthermore, for this model to succeed an unrealistically low stellar wind 
would be required. Ref. 3 suggests a phase of conservative MT from the BH progenitor 
to the companion that sets in at the end of the primary's Main Sequence, but this model 
only explains the observed masses and orbital period, failing to address the remaining 
observations. Ref. 4 proposed rotationally induced mixing as a way to keep massive stars 
from expanding significantly during the Main Sequence, preventing MT or a merger until 
the primary becomes a Wolf-Rayet star. However, this evolutionary channel increases the 
star's luminosity above that of standard models, in contrast to the observed underluminos- 
ity of the star in M33 X-7. Ref. 5 explains the observed BH's spin via a past Roche-lobe 
overflow phase from the star to the BH. However, given the extreme mass ratio between 
the components, such a phase would have evolved into a merger. 

14 
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Figure 1: Evolution of the orbital and stellar parameters of M33 X-7. From the top: masses 
(M), secondary's surface temperature [log(T cff )], secondary's luminosity [log(L 2 )], orbital 
period (P), and eccentricity (e). The different evolutionary stages are highlighted with dif- 
ferent colors: the beginning of the Main Sequence (purple), the MT phase (green), the end 
of the Main Sequence (blue), the core He burning phase for the BH progenitor (red), and 
the post BH formation phase (black) until the present time (note the non uniform x-axis). 
The grey-shaded areas represent the observational constraints as reported in Table 1 . 
The sequence comprises a ~97 M primary (Mi) and a ~32 M secondary (M 2 ) in a 
~2.9 day orbit. At the onset of the MT phase (purple/green), M x ~89.7 M , M 2 ~31 .7 M , 
and P ~3. 25 days. During the MT phase (green), the primary transfers conservatively 
~37 M to the companion (see Supplementary Information for details about the MT). 
When the system detaches (green/blue) Mi ~51 M , M 2 ~68.7M , and P ~1.6days. 
At the end of the primary's Main Sequence (blue/red) M l ~25.2 M , M 2 ~65.6M , and 
P ~2.8days. Considering evolutionary models of He stars (see Supplementary Informa- 
tion), we find that a He star of ~25M burns He at its center for ~0.38Myr and loses 
~9.1 M in its Wolf-Rayet wind. During this time the secondary loses ~1.1 M (red). Be- 
fore BH-formation (red/black) Mi ~1 6 M , M 2 ~64.5 M , and P ~3.46days. For the case 
shown, at the primary's collapse a kick of 1 20 km/s is imparted to the newly born BH (see 
Supplementary Information for the allowed kicks). 
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Figure 2: Progenitor properties and current luminosity. The circles and triangles are the 
results of detailed binary star evolution calculations for all successful sequences for a M33 
distance of 840 ± 20kpc, and of 750-1017kpc, respectively (a) possible masses for the 
progenitors; (b) possible initial orbital periods as a function of the mass ratio between the 
primary and the secondary; (c) BH X-Ray luminosity and (d) secondary's luminosity as a 
function of the secondary's mass at present. The grey and yellow shaded areas repre- 
sent the observational constraints for a distance of 840 ± 20kpc, and of 750-1017kpc, 
respectively. L x is calculated according to the Bondi and Hoyle accretion model 25 . The 
error bars are derived from the uncertainties in the stellar wind parameters (see Supple- 
mentary Information), and depict the highest and lowest L x values; they do not represent 
statistical \ a errors. The observational constraints on L 2 are calculated given the depen- 
dence of M 2 and L 2 on the distance as described in the caption of Table 1 , and accounting 
for uncertainties in the star's effective temperature, reddening, and apparent magnitude 
calculated through the ELC code. According to our model, secondaries at present more 
massive than ~65 M fail to explain the observed luminosity. Some of the data points are 
omitted for clarity. 
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This Supplementary Information provides details about the uncertainty of the distance to M33, 
the single and binary star evolution models, the binary evolution on the primary's Main Sequence 
throughout the mass transfer (MT) phase, the orbital evolution after black hole (BH) formation, 
the He stars models, the correction to the luminosity and temperature of the stellar component due 
to tidal and rotational distortions and the inclination of the system, the correction to the luminosity 
of the stellar component due to the partial-rejuvenation of the star after MT, the parameter space 
considered, and the stellar wind model used. It also contains six related figures, and additional 
references. 
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Uncertainty in the distance to M33: given the considerable uncertainty in the distance to M33 
found in the literature, we consider a compilation of 15 recent distance measurements obtained 
via different techniques (see Supplementary Figure 1). We consider the full range of uncertainty 
by combining the highest and smallest values, and obtain that the distance lies in the range 750- 
1017 kpc. 

Single and binary star evolution models: the stellar evolution models were calculated with an 
up-to-date version of Eggleton's stellar evolution code, STARS 26 ' 27,28,29 . The STARS code solves 
the equations of stellar structure and the reaction-diffusion equations governing the nuclear energy 
generation rate simultaneously on an adaptive non-Lagrangian non-Eulerian grid. We used the 
code in so-called TWIN mode, where both components of a binary are evolved simultaneously. 
This is important if both stars evolve on comparable timescale, as is the case for the M33 X-7 
progenitors. 

For the thermonuclear reaction rates we use the recommended values from NACRE 30 , with 
the exception of the 14 N(p, 7) 15 reaction, for which we use the recommended rate from Herwig 
et al. 31 and Formicola et al. 32 . The adopted opacity tables are those of Pols et al. 28 , which combine 
the OPAL opacities from Rogers & Iglesias 33 , and the low temperature molecular opacities from 
Alexander & Ferguson 34 . Conductivities come from Itoh et al. 35 and Hubbard & Lampe 36 

The assumed heavy-element composition is scaled to the solar mixture of Anders & Grevesse 37 . 



Chemical mixing due to convection ' and thermohaline mixing 9 ' is also taken into account. 
For this work we furthermore added a prescription for semi-convection based on the work of 
Langer et al. 41 All models are computed with a mixing-length to pressure scale height ratio I / H P = 
2.0. Differential rotation and mixing due to meridional circulation are not taken into account. Al- 
though mixing due to meridional circulation can be very important, the mixing efficiencies are 
also very uncertain and have to be tuned to observations. We did not wish to include extra free 
parameters in our model, which is intended to be the simplest model we can make that fits the 
observations. 

During the Main Sequence phase of the primary and for the secondary we use the mass loss 
prescription of Vink et al. 42 . When the surface H abundance by mass fraction in the primary drops 
below 0.4 we switch to the Wolf-Rayet prescription of Nugis & Lamers 43 with the metallicity 
scaling determined by Vink & de Koter 44 . If the effective temperature drops below 10,000 K the 
Vink et al. 42 prescription loses its validity and we fall back to the mass loss rate of de Jager et al. 45 . 

Binary evolution on the primary's Main Sequence, throughout the MT phase: the efficiency 
of the process of mass accretion during a phase of MT between the components of a binary system 
is still an open question, and in massive close binaries there is evidence for both quasi-conservative 
and highly non-conservative evolution 46 . Calculations of massive binary evolution have been car- 
ried out by various authors with different assumptions for the mass accretion efficiency, and the 
clear outcome of these studies is that different mass accretion efficiencies might be needed to ex- 
plain different observations. For example, Petrovic et al. 47 explored the evolutionary history of 



three of the about 20 Wolf-Rayet binary systems known in the catalogue of van der Hucht 13 , and 
they concluded that if the considered systems underwent a phase of stable MT, in order to match 
the observed systems properties a large amount of mass must have left the system. On the other 
hand, De Mink et al. 48 calculated detailed evolutionary tracks with primary masses between 3.5- 
35 M , mass ratios between the primary and the secondary components from ~1 to ~2.2, and 
orbital periods of a few days assuming both conservative and non-conservative MT. A systematic 
comparison of these evolutionary models with a sample of 50 double-lined eclipsing binaries in 
the Small Magellanic Cloud revealed that, for the 17 systems well matched by the models, "no 
single value for the efficiency of mass accretion can explain all systems", and they found good 
agreement between the model and the observed systems properties for accretion efficiencies up 
to 1 (conservative evolution). Given that we can find successful M33 X-7 formation sequences 
assuming quasi-conservative MT (with the only mass loss from the system being due to the stellar 
wind of both components) we conclude that there is no reason to invoke non-conservative MT, and 
introduce more model parameters. 

A description of some of the relevant physical quantities involved in the MT phase for the 
binary sequence described in detail in the Letter are presented in Supplementary Figure 2. 

Orbital evolution after BH formation: we study the orbital evolution of M33 X-7 from the 
moment at which the BH is formed to the present time, by examining how the orbital separa- 
tion, eccentricity, and spin frequency of the stellar component change in time. The orbital and 
spin evolution calculation accounts for the following physical mechanisms: tidal torques between 



the binary components, mass loss from the system via stellar wind, changes in the stellar radius 
during the Main Sequence lifetime of the companion star, binary angular momentum loss due to 
gravitational radiation, and accretion from the companion's stellar wind onto the BH. The BH is 
considered as a point mass and the relevant ordinary differential equations are integrated. 

The tidal torques act to synchronize the rotational motion of each star with the orbital mo- 
tion, and to circularize the orbit. The tidal evolution is calculated in the standard weak-friction 
approximation 49,50 , following the formalism of Zahn 51,49 and Hut 52 , for stars with a radiative enve- 
lope. We assume that radiative damping is the only source of dissipation. Specifically, we integrate 
numerically the set of differential equations as presented in Belczynski et al. 53 , with the only modi- 
fication being in the second-order tidal coefficient E 2 . For this coefficient, we adopt stellar models 
from Claret 54 for masses of- 63 M and ~ 79.5 M , and derive E 2 = -5.4566-7.37243-t 4 / s 0562 , 
where t MS is time in units of the main sequence lifetime. The fitting formula has only a very weak 
dependence on the initial mass for the considered range. 

Wind mass loss leads to an increase of the orbital separation and, together with the expansion 
of the star on the Main Sequence, affects the stellar spin. The evolution of the orbital separation 
driven by stellar wind is calculated following Eggleton 55 , while the change in the rotational fre- 
quency is derived assuming conservation of spin angular momentum of the star. The wind mass 
loss tends to spin up the star, while the increase in the stellar radius has the opposite effect. 

Emission of gravitational radiation acts to circularize the orbit and, together with accretion 
from the stellar wind onto the BH, shrinks the orbit. The evolution of the orbital separation and ec- 



centricity due to gravitational radiation is calculated following Junker & Schaefer 56 . The accretion 
efficiency is calculated according to Bondi & Hoyle 25 (for details about the parameters used for 
describing the stellar wind see below). For the specific case of M33 X-7, both these mechanisms 
do not significantly influence the orbital evolution. In particular, the time scale for gravitational 
radiation is longer than the timescales relative to the other physical effects mentioned above, and 
the calculated accretion efficiency is extremely small (~ 1CT 4 ). 

For each time step during the orbital evolution we calculate the Roche-lobe radius of the star 
at periastron 57 . Considering that a phase of MT via Roche-lobe overflow during this evolutionary 
stage would have been dynamically unstable, and motivated by uncertainties in the definition of 
the stellar radius inherent to different models 58 , the maximum value of the radius is set equal to the 
Roche-lobe radius. 

He star models: a common problem among most stellar evolution codes occurs when trying to 
fully remove the H envelope of a massive star that has a H burning shell right outside the He 
burning core. As the code tries to remove the last bit of the H envelope (<1% of the total stellar 
mass) it reaches inside the H burning shell, causing numerical instabilities. One way to deal with 
this problem is to dramatically increase the spatial resolution of the simulations (by a factor >10). 
However, this results in unrealistically long computational times. An alternative approach, which 
we adopt here, is to stop calculations when this phase occurs and restart from a He star model. 
This avoids the short, numerically challenging phase without losing any information. 

Using the code STARS 26,27 ' 28 ' 29 we create models of He stars with the same input physics 



used for the single and binary star evolution calculations, and with masses ranging from 3 M to 
25 M . We then evolve each model until the exhaustion of He in the core to determine the duration 
of the core He burning phase (t He ) and the corresponding amount of mass lost (AM) as a function 
of the initial mass of the He star (M He ,i). The results are shown in Supplementary Figure 3. From 
these models we derive the following relations: 



t Hc = 0.323221 + 6.24256 • M^f 762 (1) 



AM = < 



0.282 - 0.28 • M Hc ,i + 0.052 • M 2 ^ - 0.00102 • M^ e>i M He ,i < 17M 
-23.5355 + 10.262 • m(M He ,i) M He ,i > 17M 



(2) 



where t He is in Myr, and M He ,i and AM are in M . 

Correction to the luminosity and temperature due to tides, rotation and inclination: the shape 
of the star in M33 X-7 is distorted by rotation and tides. This distortion causes the temperature 
to vary over the surface of the star; the equatorial regions are colder than the poles (this is a 
consequence of the Von Zeipel theorem, which relates the effective temperature to the one-fourth 
power of the surface gravity). We use the ELC code of Orosz and Hauschildt 22 to fit a surface 
temperature map to the light curve and radial velocity curve. Using this surface temperature map, 
we find that 

(T eff ) = 0.954T polar . (3) 



Here (T cff ) is the flux-averaged effective temperature, taking into account the inclination of 74.5 
degrees reported by Orosz et al. 1 . T polar is the polar temperature; the polar temperature is the 
maximum temperature over the surface because rotation and tides have the smallest effect at the 
poles. 

The effective temperature range of 34000 to 36000 K reported by Orosz et al. 1 is a measure- 
ment of (T e ff) . However, our stellar models do not incorporate the effects of tides or rotation on the 
surface effective temperature, so a direct comparison of the model surface effective temperature 
to the measured (T cS ) is inappropriate. Because the tides and rotation have the smallest effect on 
the polar regions of the star, we choose to compare the modeled surface effective temperature with 
T po i ar = (T cff )/0.954. 

The tidal and rotational distortions of the star also have an effect on the observed luminosity. 
The true luminosity of the star is an integral over the surface of the local flux density (given by the 
Stefan-Boltzmann law): 

£truc = / oT^dA. (4) 

The luminosity quoted by Orosz et al. 1 , \og(L ohs / L Q ) = 5.72 ± 0.07, is based on the quoted visual 
magnitude, V = 18.9 ± 0.05, of the star at an inclination of % = 74.6 degrees. Given a surface 
temperature map, we can compute the luminosity that would be inferred from the orbit- averaged 
emission at an inclination of 74.6 degrees: 

A f 

L avg = - — / o-T c 4 ff cos (0 - i) dA. (5) 

-4 vis Js viB 

Here S vis is the subset of the stellar surface visible at an inclination of 74.6 degrees, A vis is the 



corresponding area, A is the total surface area of the star, and 6 is the angle of the surface normal 
with respect to the vertical axis. The cos(# — i) factor accounts for the relative orientation of the 
surface normal to the line of sight. Referring to Figure 3(b) of Orosz et al. 1 , the average visual 
magnitude of the light curve is V = 18.87; the 0.03 magnitude difference between the lightcurve 
average and the quoted visual magnitude corresponds to A log(L) = 0.01. So, we have 

log(L avg /L ) = log(L obs /L ) + 0.01. (6) 

Based on the ELC surface temperature map and equation (5), we calculate 

log(L tmc /L ) = log(L avg /L ) + 0.13, (7) 

and therefore 

log(L tme /L ) = log(L obs /L ) + 0.14. (8) 

This luminosity correction is independent of the absolute surface temperature provided the surface 
temperature profile remains fixed throughout the observed temperature range 34000 < (T eff ) < 
36000 K. The observed luminosity is lower than the true luminosity because the inclination of the 
system implies that we are looking preferentially at the colder equatorial regions of the star. When 
comparing model luminosities with observations we use L true from equation (8), not L obs from 
Orosz et al. 1 ; this amounts to a decrease in the model luminosities of 0.14. 

Correction to the luminosity due to partial-rejuvenation of the secondary: since the work of 
Hellings 59 60 on the evolution of Main Sequence mass-accreting secondaries, it is generally as- 
sumed that the accretion of matter via Roche-lobe overflow leads to so-called "rejuvenation" of 
the star. The central H abundance of the accreting secondary increases, and its internal chemical 



structure becomes almost identical to the structure of a single star of the corresponding mass. Ten 
years after Hellings, Braun and Langer 61 showed that rejuvenation does not always occur and that 
the result of mass accretion might be a star with a chemical structure unlike that of an originally 
single star. One of the most influential parameters that controls this effect is the semiconvective 
mixing efficiency that, in turn, depends on the criterion for convection used in the stellar model. 
Hellings adopted the Schwarzschild criterion, according to which the semiconvective mixing ef- 
ficiency is infinite, while Braun and Langer used the Ledoux criterion with a finite value for this 
parameter. Braun and Langer showed that, despite the increase in luminosity as a result of mass 
accretion, the non-rejuvenated models appear to be underluminous for their new mass during the 
remaining Main Sequence evolution. Following Braun and Langer, we use the Ledoux criterion for 
convection in our detailed single and binary evolution calculations, using primordial composition 
for the transferred material, and we adopt for the semiconvective efficiency parameter a value of 
a sc = 0.0025 (which was calibrated using some of the results reported by Braun and Langer 61 ). 
Our results confirm that the rejuvenation of the secondary component after MT was at most partial 
(see Supplementary Figure 4 for an example). 

The parameter space: keeping in mind that the short observed orbital period cannot be the re- 
sult of a common envelope phase from the BH progenitor to the companion star, we explore the 
evolution of binary systems that start their life already in a tight orbit, hence undergoing a MT 
phase during the core H burning phase of the primary component. We perform the binary evolu- 
tion calculations from the Zero-Age Main Sequence until the end of the primary's Main Sequence 
considering various combinations of initial masses and orbital periods. Specifically, we evolve pri- 



maries and secondaries with masses between 20-130 M and 10-100 M , respectively, and initial 
orbital periods ranging from 1 to 10 days. Based on the observed masses of the two components 
we use a different density of models in different regions of the parameter space. 

Since we consider MT during the core H burning phase of the primary, we perform a first 
scan of the data rejecting the sequences where the primary overfills its Roche-lobe only after the 
end of its Main Sequence. Furthermore, given the high mass of the BH companion, we also 
exclude the sequences where the secondary ends up transferring mass to the primary after having 
accreted mass from it. After the end of the primary's Main Sequence, given that the star is a He 
star we use equations (1) and (2) to calculate the amount of mass lost from the binary components, 
and the corresponding change in orbital period during the core He burning phase of the primary 
until collapse. We reject the sequences where the masses of the two stars at the end of the core 
He burning phase of the primary are lower than the minimum observed values for the minimum 
distance to M33 of 750 kpc. Since no episodes of Roche-lobe overflow can have occurred from the 
collapse of the primary until the present time, we then evolve each secondary as a single star, and 
exclude the sequences where the mass of the star does not fall within the observed range when the 
model matches the observed effective temperature and luminosity. After the primary's collapse, we 
consider a variety of orbital configurations by scanning the parameter space made up of the kick 
magnitude (14), orbital separation (a pos tBH), and eccentricity (e post BH)- Specifically, we consider 
isotropic kicks between 0-1300 km/s, orbital separations between 0-100 R , and eccentricities 
between and 1. The requirement that the system must remain bound after the BH is formed, and 
that the direction of the kick must be real impose constraints on the pre- and post- BH formation 



orbital parameters 62 63 . We then study the orbital evolution after the formation of the BH, and 
interrupt the calculation when the orbital period crosses the observed value; then the eccentricity 
of the orbit must fall within the observed la range. 

Finally, of the sequences that fulfill all the above requirements, we compare the stellar radius 
at present with the distance from the center of the star to the point through which mass would 
flow from the star to the BH in case of Roche-lobe overflow 64 . We reject the sequences where the 
secondary companion is transferring mass to the BH at the present time. Supplementary Figure 5 
shows the masses of the components at present according to our model. Supplementary Figure 6 
shows the allowed BH kicks, BH progenitor masses, orbital separations and eccentricities post-BH 
formation, for all the successful sequences. 

According to our model, the BH progenitor mass lies within the observed range for the BH 
mass, and, hence, no baryonic mass is ejected at collapse. Furthermore, the allowed eccentricities 
post-BH formation are constrained to be between 0.012-0.026. On one hand, the lack of mass 
ejection at BH formation, and the small induced eccentricity, could imply that the BH did not 
received a high kick at formation. On the other hand, due to the lack of kinematic information, we 
can not exclude kicks as high as ~850km/s. This apparent discrepancy is explained by the fact 
that the kick is constrained to point mostly orthogonal to the orbital plane. In this case, a higher 
kick results in a more tilted orbit, while the orbital eccentricity does not change significantly, but 
enough to explain the observed eccentricity of 0.0185±0.0077. An upper limit to the magnitude 
of the kick is given by the observationally inferred positive spin of the BH. Kicks higher than 



~850km/s would flip the orbital plane and that would result in a negative value for the BH's spin. 



Stellar wind model and X-Ray luminosity: to determine the stellar wind properties that enter the 
Bondi & Hoyle 25 accretion model we follow Lamers & Cassinelli 65 , and we adopt a velocity law 
of the type 

v (r) = v + (v oo -v )(l--) , (9) 



r j 

where v is the escape velocity at the stellar surface, is the velocity of the wind at infinity, 
R is the photospheric radius, and f3 is an index which typically ranges from 0.8 to 1.2 (we use 
f3 — 1.0 ± 0.2). The escape velocity is defined as 

"2GM(i -ry 1/2 



v 



R 



(10) 



where 



r = 7.«xio-v.(A)(£). ( „) 

For the electron scattering coefficient a e we use 

a c = 0.40l|±f , (12) 
1 + 3e 

where q is the fraction of He ++ , and e = He/ (H + He). Following Lamers & Leitherer 66 , we 
adopt e = 0.15 ± 0.05, which is appropriate for an O-type stars of spectral class III, q = 1 if the 
effective temperature of the star (T cff ) is > 35, 000 K, or q = 1/2 if 30, 000 < T cff < 35, 000. For 
the velocity of the wind at infinity we adopt Voo/vq = 3.085 ± 1.075, which is again appropriate 
for a star of the spectral class observed for the companion star in M33 X-7 67 . To calculate the 
mass accretion rate via stellar wind according to the Bondi & Hoyle accretion model, we follow 



Belczynski et al. 53 (and references therein): 



\ 2 r 

^wind ( GM &CC \ Q^wind -^don.wind 



M acc , wind - ^ y ^ j ^ 2 (i + y2)3/2 . (13) 

We use a wind = 3/2, and F wind = 1. If M acCiWind exceeds 0.8M doiljWind , F wind is set such that 
Ma,cc,wind = 0.8M don ,wind- Mice is the mass of the accreting component, M doil]Wind is the mass 
loss rate via wind, and V 2 = V 2 cc olh /V 2 ind . The orbital velocity of the accretor is given by 
Kcc,orb = C(M acc + M don ) / a, and V^ ind is the squared of the velocity of the wind at r = a, where 
a is the orbital separation. Following Belczynski et al. 53 , the bolometric luminosity is calculated 
from the mass accretion rate 

boi = e , (14) 

-'•-ace 

where e gives a conversion efficiency of gravitational binding energy to radiation associated with 
accretion onto a compact object, and is equal to 0.5 for accretion onto a BH. R acc is the radius of 
the accretor, and we calculate it from the observationally inferred spin following Bardeen et al. 68 . 
Finally, the X-Ray luminosity is calculated from the bolometric accretion luminosity via 53 : 

Lx = VboiLhoh (15) 

where we use rj = 0.8 ± 0.1. 
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4 Supplementary Figures And Legends 
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Supplementary Figure 1 Recent measurements of the distance to M33. From the 
bottom, 1 : U et al. 6 , 2: Scowcroft et al. 69 , 3: Orosz et al.\ 4: Bonanos et al. 70 , 5: Sara- 
jedini et al. 71 , 6: Ciardullo et al. 72 , 7: Galleti et al. 73 , 8: McConnachie et al. 74 , 9: Tiede 
et al. 75 , 10 and 11: Kim et al. 76 , 12: Lee et al. 77 , 13: Freedman et al. 78 , 14: Pierce et 



al. 79 , 1 5: Sarajedini et al. 80 . The last 1 2 estimates have been calculated from the distance 
moduli listed in the compilation given by Bonanos et al. 70 . We do not include the distance 
measured from Brunthaler et al. 81 because of its big uncertainty. In fact, both the sta- 
tistical and systematic errors are quite large compared to other methods. A longer time 
baseline of observations would help reduce the statistical error (Bonanos 2010, private 
communication). 




Supplementary Figure 2 The conservative MT phase. From the top: masses of the 
components, mass rate of change (M), typical stellar evolution timescales {TS), Roche- 
lobe overflow filling factor (R/R L ), stellar radius (R), and primary's H surface abundance 
by mass fraction (X s ) as a function of time for the evolutionary sequence described in the 
Letter. Black and red lines indicate the primary and secondary component, respectively 



M tr is the mass transfer rate. The MT (or Roche-lobe overflow) phase is denoted by a 
positive value for log (R/R L ), where R L is the star's Roche-lobe radius (the dotted line 
in the corresponding plot is given as a reference for log (R/Rl) = 0). When X s drops 
below 0.4 the Wolf-Rayet regime is entered (the dotted lines is given as a reference to this 
value). Given the short duration of the MT phase with respect to the total Main Sequence 
lifetime of the primary star, we use a non uniform x-axis. The ~97 M primary evolves 
and expands faster than the secondary, and after ~1.8Myr overfills its Roche-lobe and 
begins transferring part of its envelope to the companion. Within the first ~37,000yr into 
the MT phase the orbital period decreases and, as a result, the rate at which the primary 
transfers mass increases from ~10 8 M /yr to ~3- 10 3 M /yr. At this time, the MT and 
mass accretion timescales (M 1 /M tT and M 2 /M tr , respectively) become comparable to 
the stars' thermal timescales. This brings the components out of thermal equilibrium (the 
components remain in hydrostatic equilibrium). The secondary overfills its Roche-lobe 
radius as well for a short time (~9,000yr) without transferring any mass, but when the 
component masses become equal and the orbit begins expanding, both stars recover their 
thermal equilibrium, and the secondary detaches. Once the secondary is detached and 
the orbital period (and Roche-lobe radius) are increasing, the primary keeps transferring 
mass for ~59,000yr at a decreasing rate. When X s of the primary drops below 0.4 the 
star enters the Wolf-Rayet regime, and the corresponding stronger stellar wind interrupts 
the MT, and removes the remaining stellar envelope to expose the He core. 
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Supplementary Figure 3 Core He burning phase for He stars with different initial 
masses. Top: duration of the core He burning phase as a function of the initial mass. 
Bottom: amount of mass lost during the core He burning phase as a function of the initial 
mass. The dots represent He star models, while the solid lines are the fits in equations 
(1) and (2). 




Supplementary Figure 4 Non-rejuvenation. Hertzsprung-Russell diagram for a 32 M 
star that accreted 37 M , and for a single star that starts its life with 69 M . The stellar 
model that underwent mass accretion is the secondary component described in detail 
in the Letter. The luminosity does not include the correction due to tidal and rotational 
distortions and the inclination of the system with respect to the line of sight. 




Supplementary Figure 5 Masses of the components at present. The circles and 
triangles are the results of detailed binary star evolution calculations for all successful 
sequences for a distance to M33 of 840 ± 20kpc, and of 750-1017kpc, respectively. 
The grey and yellow shaded areas represent the observational constraints for a distance 
of 840 ± 20kpc, and of 750-1017kpc, respectively The observational constraints are 
calculated given the dependence between the masses of the components M BH = 6.19 + 
0.13 -M 2 , and accounting for uncertainties in the star's effective temperature, reddening, 



and apparent magnitude calculated through the ELC code. Some of the data points are 
omitted for clarity. 
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Supplementary Figure 6 Pre- and post- BH formation orbital parameters. The cir- 
cles and triangles are the results of detailed binary star evolution calculations for all suc- 
cessful sequences for a distance to M33 of 840 ± 20kpc, and of 750-1017kpc, respec- 



tively. (a) orbital eccentricity as a function of the orbital separation post-BH formation; 
(b) mass of the BH progenitor as a function of the kick magnitude; (c) change in the or- 
bital inclination at BH formation as a function of the kick magnitude. The BH progenitor 
mass accounts for the 10% of rest mass energy that is released as the BH's gravita- 
tional energy at collapse. According to our model, the BH progenitor mass lies within 
the observed range for the BH mass (15.65 ± 1.45M for a distance of 840 ± 20kpc, 
and between 1 3.5-20 M for a distance of 750-1017kpc). Our model allows kicks from 
~ 10 km/s to ~ 850 km/s. Some of the data points are omitted for clarity. 



